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ABSTRACT 


Model validation and flight test data analysis require careful consideration of the effects of uncertainty, 
noise, and nonlinearity. Uncertainty prevails in the data analysis techniques and results in a composite 
model uncertainty from unmodeled dynamics, assumptions and mechanics of the estimation procedures, 
noise, and nonlinearity. A fundamental requirement for reliable and robust model development is an at- 
tempt to account for each of these sources of error, in particular, for model validation, robust stability 
prediction, and flight control system development. This paper is concerned with data processing proce- 
dures for uncertainty reduction in model validation for stability estimation and nonlinear identification. 

F/A-18 Active Aeroelastic Wing (AAW) aircraft data is used to demonstrate signal representation 
effects on uncertain model development, stability estimation, and nonlinear identification. Data is decom- 
posed using adaptive orthonormal best-basis and wavelet-basis signal decompositions for signal denoising 
into linear and nonlinear identification algorithms. Nonlinear identification from a wavelet-based Volterra 
kernel procedure is used to extract nonlinear dynamics from aeroelastic responses, and to assist model 
development and uncertainty reduction for model validation and stability prediction by removing a class 
of nonlinearity from the uncertainty. 


INTRODUCTION 


Aeroelastic and aeroservoelastic (ASE) flight data analysis and model validation must necessarily cope 
with uncertainty, noise, and nonlinearity. Model uncertainty between a high-fidelity computational model, 
reduced-order (ref. [1]) model, or a linear model, and a descriptive representation of the actual aircraft 
dynamics is inevitable. ASE systems comprise interactions of generally multi-input multi-output sampled- 
data control feedback with actuation dynamics coupled with ASE. Highly augmented closed-loop flight 
test data require extra care in distinguishing system component dynamics. Important considerations of 
uncertainty classification relative to the signal processing and identification procedures, versus an actual 
model mis-match, must be assessed and minimized for accurate model updates. 

Identification of ASE dynamics presents one of the greatest challenges for flight data analysis. Un- 
steady aerodynamics, turbulence, shock effects, and spatially local behaviors are often difficult to quan- 
tify, especially with inadequate instrumentation. Exogenous inputs and asymmetric characteristics create 
problems for observability and identifiability. Model verification over an extensive flight envelope presents 
more challenges. Test data acquisition is expensive so maneuvers are designed for maximum efficiency and 
data quality with sometimes costly consequences. Verification methods are desired which accurately and 
efficiently include identification of critical parameters, address mismodeling and unmodeled dynamics, 
deal with test condition and system variability, and derive data-consistent parametric and nonparametric 
uncertainty descriptions (ref. [2, 3]). Non-parametric descriptions are addressed in this paper. 

This paper attempts to approach estimation of ASE models with procedures to account for uncertainty, 
noise, and nonlinearity. Set membership identification in the form of bounded error estimation (ref. [2, 4]) 
has been used previously to characterize feasible sets of parameters with uncertainty estimates consistent 
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with the data, model structure, and prior information on uncertainty bounds. Minimum upper error bounds 
computed with parameter estimates result in a feasible set described as a function of the error bounds, 
thereby resulting in error bounds describing a smallest non-empty feasible set. These error bounds can 
be refined with appropriate data decompositions which extract relevant features from a time-frequency 
perspective. 

Extraction of local time and frequency content of signals has been shown to be useful in nonstationary 
signal analysis for classification, recognition, and interpretation (ref. [5, 6]). Automated techniques are 
being pursued to detect impulsive, wide-band, and harmonic components by adaptive tracking of time- 
frequency and time-scale variations of multi-component signals (ref. [7, 8, 9]). Best basis algorithms seek 
to generate fast adaptive sparse representations that have enhanced resolution compared to non-adaptive 
approaches. These atomic decompositions act on overcomplete dictionaries for more general analysis over 
a wide class of signals (ref. [10, 11]). However, adaptive time-frequency localization is more efficient with 
orthogonal or approximately orthogonal data decompositions that project signals onto functions that form 
an orthonormal basis of finite energy signals (ref. [12, 7, 13]). 

Application of an adaptive orthogonal data decomposition is used here to optimally extract informa- 
tion content. One technique uses a fast adaptive best-basis expansion of signals in local trigonometric 
bases (ref. [12, 14])(LTBs). Another near-orthogonal redundant wavelet basis (ref. [7, 15, 16, 17]) is also 
used for nonstationary time-frequency representations. Transfer functions are derived to demonstrate a 
less conservative uncertainty bound for robust stability analysis when using these signal representations. 

These signal representations are incorporated in a nonlinear multiwavelet-based Volterra kernel iden- 
tification procedure (ref. [18, 19]) to discern nonlinear-plus -noise dynamics from a best underlying lin- 
ear representation. Kernels are identified from wavelet-processed data and discussion is presented about 
the applicability and feasibility of identifying Volterra kernels from nonlinear F/A-18 Active Aeroelastic 
Wing (ref. [20]) flight data. 


AEROSERVOELASTIC DATA ANALYSIS 


The ASE open-loop plant model includes rigid body and elastic modes, coupled high-order actuator 
dynamics, and control surface modal dynamics (ref. [3]). Including the aerodynamic lag states, the ASE 
state equations are 


x = Ax + Bu 1 y = Cx + Du ; x = \r]sr) r r] e r] r r) e r) a 5 J 

consisting of input control surface commands u, actuator states rjs, rigid body states rj r , flexible mode 
states Tj e , aerodynamic lag states r/ a , and control surface displacements S. Aeroservoelastic plant, P, 
is therefore represented as the state-space operator. Associated with this time-domain representation is 
the transfer function, P(s), a function of the complex Laplace variable, s, such that y = P(s)u, where 
P(s ) = D + C(sl — A)~ l B. Controller K(s) is modeled similarly, but being a digital implementation, it 
is modeled as a function of discrete complex variable, z, as K(z = e sT ) specified by the sampling time T 
and a zero-order sample-hold at the input of the controller. 
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Figure 1 : Closed-loop ASE model with additive and multiplicative plant uncertainty. 


Multiplicative uncertainty is used to represent unmodeled dynamics and errors from the feedback sen- 
sors ( w 0 = A 0 z 0 ) and actuator input commands (vjj = as shown in Figure 1. Each of A 0 and A* 

are diagonal complex perturbations of appropriate output or input dimensions, with suitable weights W t 
such that the uncertainty pertubations are norm-bounded by unity. Multiplicative perturbation at the output 
results in perturbed model P 0 = (I + A 0 )P, and at the input the perturbed model is Pi = P(I + A*). 
Structured robust stability tests representing these types of uncertainty are described from the comple- 
mentary sensitivity matrix functions of complex P(s) and K(s), namely, T 0 = PK(I + PK)~ l and 
Ti = I\P(I + KP)~ l . Another uncertainty characterization is the additive perturbation, A a , with per- 
turbed plant P a = P + A a , which is depicted in Figure 1 (w a = A a z a ). The effect of real parametric 
uncertainty at the plant input or output can be shown to represent an additive uncertainty in the plant trans- 
fer function (ref. [2]). Necessary and sufficient conditions for structured robust stability are then derived 
from the following tests on loop sensitivity functions (ref. [21, 22]) (where HA^oo < 1, 1 1 A* 1 1 oc < 1, and 
1 1 A a | |oo < 1) with the structured singular value, ju, depending on the type of uncertainty, {A*, A 0 , A a }. 

Ma 0 (W 2 T 0 W 1 ) < 1, <UAi(W 4 TiW 3 ) < L <u Aa (W 6 K(I + PK^Ws) < 1. 


Data Decompositions 


An uncertainty bound estimation technique based on robust minimax estimation (ref. [2, 4]) is used 
here to demonstrate the effects of signal representations on robust stability analysis. Minimum upper error 
bounds are computed with the ASE model parameter estimates such that the feasible estimation set is 
described as a function of the error bounds with no apriori assumption on the noise or model uncertainty. 
Transfer functions are computed in the robust minimax estimation scheme to identify input-output param- 
eters and structured error bounds consistent with the data. These uncertainty estimates provide structured 
uncertainty representations for robust stability analysis useful for model validation and control system 
design. 
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Two types of data decompositions (projections of the data onto a set of pre-defined basis functions) 
to be discussed in this paper are of concern for either linear or nonlinear signal analysis and system 
identification. Many other options (ref. [10, 8, 9, 11, 23]) based on time-frequency filtering, localization, 
and wavelet denoising (ref. [13]) ideas have shown promise for signal enhancement, interpretation, and 
analysis, and will be pursued for ASE flight test dynamics in future research. The two chosen in this paper, 
(1) local trigonometric (ref. [12]) and, (2) Morlet wavelet (ref. [15, 16]) bases, are illustrative of other 
methods and have sound theoretical and practical justification. 

Local time-frequency signal decompositions generally require basis functions well localized in both 
time and frequency. Orthogonality helps for algorithmic efficiency. An adaptive decomposition is desired 
for better accuracy and sparsity in the representation. Signals are projected onto orthonormal bases each 
with different time-frequency properties before choosing the best basis with some (information cost) op- 
timality criteria. Popular and fast algorithms for best-basis selection include wavelet packets (ref. [14]) 
which allow for a general basis representaion but usually require a dyadic decimation at each scale, and 
it can be difficult to relate coefficients between scales or from a scale to the original signal. LTBs have 
excellent time and frequency resolution and allow for simple and fast algorithms (ref. [12, 13, 14]) while 
relying only on trigonometric functions, which may not be appropriate for signals of interest. 

LTBs assign a ramp, or rising cutoff function r(t) = exp[jp(t)] sin[0(t)], which rises from zero to 
unity in the interval {—1, 1} for p G {2n7r, 2m7r} and 9 G {0, |}. Windows are created from the ramp 
functions applied to both sides of each interval. Local cosines are constructed by dilating, normalizing, and 
translating windowed block cosines over a disjoint partition of intervals of the signal. Cosines and/or sines 
are then used as bases to orthogonally project the data onto compactly supported windowed trigonometric 
functions which are well-localized in time and frequency. Adaptivity depends on best-basis algorithms 
choosing optimal partitions and corresponding coefficients for the expansion. A reasonable criteria is 
the decay rate of coefficients as a function of the partition, which is quantified by the derivatives locally 
at the intervals. Based on this criteria the decision to divide intervals is made locally independent of 
other changes (ref. [12]). By using the Lourier transform property of derivation, < >= 

(— j7rn/|/|) m < /, e _J7rn /l 7 l > for partition length \I\, the algorithm is efficient. 

Wavelet decompositions also have desirable properties for signal decompositons if the wavelet basis 
functions have properties inherent to the signal. Morlet wavelet multiresolution analysis can provide 
efficient (well-matched to the signal) and descriptive coefficients for feature detection and identification 
purposes (ref. [7, 24, 16]) much like the LTBs for sinusoidal-like signals. 


Results 


Llutter and ASE clearance in both the subsonic and supersonic envelopes for the Active Aeroelastic 
Wing is accomplished with multisine sweeps summed to pertinent actuator commands, as was done previ- 
ously on the NASA Dryden L/A-18 Systems Research Aircraft (ref. [2]). The 56 structural accelerometers 
are distributed on the fuselage, empennages, and all control surfaces. Control system feedbacks, com- 
mands, and surface positions are also available for ASE monitoring. Demonstration of the effect of signal 
decompositions into the LTBs and the Morlet wavelet functions will now be demonstrated. 
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This paper is concerned with ASE stability estimation and aeroelastic dynamics identification with 
signals represented as in Figure 2. Here is the right rear-wingtip accelerometer response at MachO. 85, 
15,000ft altitude during a multisine input to the collective aileron command. The sweeps run from 3- 
35Hz for 26sec, and this segment is a portion between 5-10sec during excitation of the first symmetric 
wing bending mode. The top plot is the raw data at 200sps and the second plot is LTB -filtered with a 
hard denoising threshold of the LTB coefficients at 90% (ref. [12]). Hard thresholding corresponds to an 
optimization of the signal into a best basis by minimizing the residual norm amongst all data sets with 
a fixed number of nonzero coefficients determined by the threshold (ref. [11]). Despite this extreme 
threshold, the signal looks remarkably similar to the original. However, the LTB -filtered signal must 
correspond to a trigonometric basis, and the residual shown on the same plot (second from the top) is 
white noise. 
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Figure 2: AAW right rear-wingtip accelerometer response (g’s) of first symmetric wing bending mode from multi- 
sine excitation to the symmetric aileron command. From top to bottom: 

Raw data, LTB-filtered with residual, Morlet-filtered, and residual of Morlet-filtered signal. 

The Morlet-filtered data in the third plot of Figure 2 is clearly a sinusoidal basis of the data in-phase 
with the original signal. A residual of the Morlet wavelet coefficients in the fourth plot is calculated as the 
signal reconstruction from the coefficients that are not used in the reconstruction of the filtered signal of 
the third plot. These residual wavelet coefficients are extracted from the scalogram (time-frequency map 
of the wavelet coefficients) outside the region of the input signal. Hence, these contributions are due to 
extraneous signal content outside of the commanded input multisine sweep such as from exogenous input, 
nonlinearity, and noise in the response. These Morlet-derived signal representations will assist in linear vs 
nonlinear dynamics identification using Volterra kernels in the next section. 
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Figure 3: MachO. 85, 5,000 feet; Additive perturbation analysis (MaxAdd), complex multiplicative (Complex), and 
nominal (Nominal) model analysis using complentary sensitivity functions T t and T 0 . 

Left plots are longitudinal analyses, and right plots arc lateral-directional analyses. 

Raw data (top), LTB-filtered data (middle), and Morlet-liltered data (bottom). 


To demonstrate the effects of these signal representations on stability analysis, the raw data, LTB- 
filtered data, and Morlet-filtered feedback responses and control commands are used to compute the com- 
plementary sensitivity functions {T*, T 0 } and perform the structured robust stability analyses represented 
by Figure 1 and expressions (1). In the top plots of Figure 3 and Figure 4, raw data from two flight con- 
ditions are chosen for longitudinal and lateral-directional robust stability analyses. These conditions are 
chosen as a transonic (MachO. 95, 15K feet) and a highly turbulent (Mach0.85, 5K feet) case. Multisine 
inputs are aileron, outboard leading edge, collective stabilator, and rudder commands, and outputs are nor- 
mal and lateral acceleration, roll rate, pitch rate, and yaw rate feedbacks. Additive perturbation analysis 
(MaxAdd) uses A a from Figure 1 derived from the data, complex multiplicative (Complex) analysis is de- 
rived from updated transfer functions and {Aj, A 0 } from (1) based on corrections to the nominal models 
to match the data in a minmiax optimization procedure (ref. [2, 4]), and the nominal (Nominal) analysis is 
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Longitudinal, MachO. 95 Lateral-directional, MachO. 95 



Frequency, Hz Frequency, Hz 

Figure 4: MachO.95, 15,000 feet; Additive perturbation analysis (MaxAdd), complex multiplicative (Complex), 
and nominal (Nominal) model analysis using complentary sensitivity functions T t and T 0 . 

Left plots are longitudinal analyses, and right plots arc lateral-directional analyses. 

Raw data (top), LTB-filtered data (middle), and Morlet-filtered data (bottom). 

also complex multiplicative uncertainty analysis at the input and output from the baseline analysis models. 

Note the large errors between nominal and updated models at MachO. 85, and high gain levels in the 
lateral-directional responses near 10Hz (near and above OdB, and large uncertainty). Since the time- 
domain data does not indicate a tendency towards resonance or instability, the responses are suspect due to 
turbulence levels or transonic effects. LTB-filtered data is used to generate the responses in the middle plots 
of Figures 3 and 4, where there is much improvement in the longitudinal uncertainty levels and matches 
between updated and nominal models at MachO. 85, and some improvement in the still-conservative gain 
levels of the lateral-directional responses. 

Finally, in the bottom plots of Figure 3 and Figure 4 it is shown that the Morlet-filtered results gener- 
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ally produce less conservative stability margin and more realistic lateral-directional gain levels where the 
modal peaks near 10Hz are more reasonable with reference to the flight data. This shows a reasonable im- 
provement in the Morlet-filtered responses especially in the transonic and lateral-directional comparisons 
with raw and LTB- filtered data. 


NONLINEAR AEROELASTIC DATA ANALYSIS 


Volterra series representations provide a convenient framework for the analysis of nonlinear dynamical 
systems. The Volterra theory of nonlinear systems states that the system output y can be expressed in terms 
of an infinite series of integral operators of increasing order (ref. [25]). In practice, the series is truncated 
and this paper considers Volterra models that include only the first and second-order operators. For a 
causal, time-invariant, single-input/single-output system, the first and second-order Volterra operators take 
the form 

Vi(t)= [ hi(£)u(t - £)d£ (1) 

Jo 

Vi (t) = f [ h 2 (i,r])u(t - £)u(t - rj)d£dri (2) 

Jo Jo 

where u is the input and hi, h -2 are the first and second-order Volterra kernels. Collectively, the Volterra 
kernels provide a model of the system since, once the kernels have been identified, the response to any 
arbitrary input can be determined. The first-order kernel represents the linear dynamics of the system while 
the kernels of second-order and higher characterize the nonlinear dynamics. It should be noted that, for a 
linear system, the first-order kernel is equivalent to the impulse response of the system and the output is 
given by equation (1). Therefore, the Volterra theory can be viewed as an extension of the concept of linear 
convolution to nonlinear systems. In general, Volterra series are applicable to systems that exhibit fading 
memory. Fading memory asserts that past inputs have a diminishing influence on the present output. This 
implies that all of the Volterra kernels of a given system decay to zero in a finite period of time. 


Volterra kernel identification 


The difficulty in using Volterra series to model dynamical systems lies in the identification of the 
Volterra kernels. Fundamentally, kernel identification is an ill-posed, inverse problem since the objective 
is to determine the system model from input and output measurements. Also, a large number of coefficients 
are often needed to represent Volterra kernels, with the number increasing geometrically with the order of 
the kernel. For example, a discrete Volterra model requires on the order of N p coefficients to represent the 
pth-order kernel for an ^-dimensional data set. In view of these limitations, many different approaches 
have been taken to identify Volterra kernels in both the time and frequency domains. Direct measurement 
techniques have been employed whereby the kernels are identified in terms of the response of the system 
to specific inputs. These include harmonic probing and the application of impulsive inputs to discrete-time 
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systems (ref. [26]). Many approaches are statistical in nature and can be viewed as variations of the cross- 
correlation technique developed by Lee and Schetzen (ref. [27]) for Gaussian white noise inputs. Other 
methods, including the approach taken in this paper, seek to expand the kernels in terms of a relatively 
small number of basis functions such as Laguerre polynomials (ref. [28]) or decaying exponentials (ref. 
[29]). 

In this paper, orthonormal, piecewise-polynomial multiwavelets are used to approximate the first and 
second-order kernels of nonlinear aeroelastic systems. Wavelets are compactly- supported, oscillatory 
functions that are constructed to satisfy certain properties such as orthogonality, smoothness, and sym- 
metry requirements (ref. [13]). Multiwavelets compose a set of wavelet functions {if) 1 , . . . , ip r } that are 
generated from a set of scaling functions {(f) 1 , . . . , cf) r } (ref. [30]). The scaled translates and dilates of 
the multiwavelets form a basis for L 2 ( R), the space of square-integrable functions. The scaling functions 
generate a multiresolution analysis, a series of nested approximation spaces {Vj}j e z- Each space Vj is 
defined as 

Vj := span {(V ]k : k G Z, s G {1, . . . , r}} 

where 

- k) 

Hence, increasing the dilation index j by one effectively doubles the resolution of the space Vj. The 
translation index k gives the position of the function on the real line. The multiwavelets are constructed 
so that they span the orthogonal complement spaces {Wj}: 

Wj = Vj + 1 © Vj, j G Z (3) 

Each wavelet space Wj is defined as 

Wj := span : k G Z, s G {1, . . . , r}} 

where 

Vi* = 2"V(2’ - k) 

Applying equation (3) recursively, the approximation space Vj +1 can be decomposed as 

Vj . |-i = Wj © Wj — i © • • • © W j0 © Vj 0 (4) 

where jo is the coarsest resolution level in the decomposition. 

This paper employs piecewise-quadratic multiwavelets that have been constructed using the technique 
of intertwining (ref. [30]). This process derives four scaling functions and associated multiwavelets from 
the classical quadratic finite element basis functions. The details of this construction and plots of the 
multiwavelets are given in Ref. [18]. This class of multiwavelets is well-suited for the approximation of 
Volterra kernels because the functions are orthogonal, symmetric or antisymmetric, and are easily adapted 
to the finite domains over which the kernels are supported. 

In using these multiwavelets for kernel identification, the first-order kernel is first written in terms of a 
fine-resolution approximation space Vj as 

4 

mo = EE<A® < 5 > 

s=l k 
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where {o;^} are constant coefficients. Equation (4) implies that an equivalent, multiscale representation 
of the kernel can be obtained as 


Mf) = EEE + E E 

l=jo s=l k s=l k 

This multiscale representation, termed the discrete wavelet transform, is desirable because it decomposes 
the kernel into a coarse-scale average over V 3Q and a series of details of varying resolution over the wavelet 
spaces Wj o: . . . , Wj_ i. Furthermore, many of the coefficients in the multiscale expansion are often close 
to zero and can be neglected, leading to a sparse representation. The kernel identification problem is easier 
to formulate in terms of the single-scale representation in equation (5). Fast, recursive algorithms are then 
applied to relate the single-scale coefficients to the multiscale wavelet coefficients. 


In a similar manner, the second-order kernel is expanded in terms of two-dimensional scaling func- 
tions and multi wavelets. These two dimensional functions are obtained as the tensor products of the 
one-dimensional scaling functions and wavelets. A single-scale approximation of the second-order kernel 
can be written as 

4 


= I] ^ a Uk,rn)^AOrt, m (.V) 


s,i— 1 k,m 


(6) 


and the equivalent multiscale expansion is given by 


^) = EE wiv) + 

s,i— 1 k,m l— jo s,t=l k,m 

+ 

l— jo 5,t— 1 k,m l— jo s,t= 1 k,m 


(7) 


It has been shown that the Volterra kernels can be assumed to be symmetric. Under this assumption, the 
number of unique coefficients in equations (6) and (7) can be reduced by almost a factor of two. 


The multiwavelet representations of the first and second-order kernels are substituted into equations 
(1) and (2), respectively. Then, using a zero-order hold approximation of the input and output, the kernel 
identification problem reduces to the matrix equation 

(see Ref. [18] for details). In equation (8), y. represents a vector of discrete outputs that have been 
sampled at 2 J sps and y i , y 2 . are vectors of discrete first and second-order outputs. The vectors and 
are composed of the multiscale wavelet coefficients that represent the first and second-order kernels. 
Equation (8) is solved, in a least-squares sense, for the first and second-order kernel coefficients. In many 
cases, the vectors /3 and ft 2 can be truncated to obtain reduced-order representations of the kernels. 
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Results 


The multiwavelet-based Volterra kernel identification method was used to extract first and second- 
order kernels from AAW flight data. In previous work, a similar approach was applied to flight data 
from the Aero structures Test Wing (ATW), a small-scale test wing that was designed to study flutter (ref. 
[19]). In that study, first-order kernels were successfully identified from the flight data at various flight 
conditions. The identified second-order kernels were essentially zero, indicating that the ATW was a 
predominantly linear system with two degrees-of-freedom. The AAW, on the other hand, is a much more 
complex system with many coupled degrees-of-freedom, and the flight data clearly indicates the presence 
of aeroelastic nonlinearities. 

In practice, first and second-order Volterra kernels must be identified simultaneously because the terms 
in the Volterra series are not orthogonal. The Morlet wavelet filtering procedure can be used to simplify 
the identification procedure by decomposing the data into a linear component and a residual. The linear 
component is generated by correlating the input and output frequencies in the time domain. The portion 
of the output that matches the input frequencies is assumed to be linear. The residual data is assumed to be 
composed of nonlinear data and noise. In this manner, a good approximation of the linear and nonlinear 
portions of the response data is obtained. As an additional benefit, much of the noise in the data is removed 
during the Morlet filtering process. The first-order kernel is then identified from the Morlet-filtered data, 
and the second-order kernel is identified from the residual data. Of course, it is certainly possible that the 
Morlet-filtered data contains some nonlinear response at the input frequencies. Such response could be 
modeled in terms of odd-order kernels of third order and higher, which are not considered in this paper. 

As an example, first and second-order kernels were extracted from AAW flight data at a flight condition 
of 15, 000 ft, Mach number 0.85. The input was a multisine collective aileron sweep and the output was 
taken as accelerometer data from the forward right wing just inside the wing fold. This data, which was 
sampled at 128 sps, is shown in Figure 5. It should be noted that the mean of both data sets has been 
removed to ignore variations in the trim condition during the maneuver. The Morlet filtering procedure 
was applied to the accelerometer data. The Morlet-filtered response and the residual are shown in Figure 
6 . 


A first-order kernel was identified from the Morlet-filtered data, which was assumed to be linear. 
Then, a symmetric, second-order kernel was extracted from the residual data, which was assumed to be 
composed of nonlinear data and noise. The identified kernels are shown in Figure 7. The memory, or time 
duration, of each kernel was chosen by trial and error to be the minimum time for each kernel to decay 
to zero. Minimizing the kernel lengths reduces the number of coefficients that must be identified. This is 
especially significant for the second-order kernel, where doubling the memory would increase the number 
of coefficients by a factor of four. Similarly, the number of wavelet coefficients needed to represent each 
kernel, which is directly related to the number of levels of resolution used, was also determined via trial 
and error. Basically, the number of resolution levels was increased until the kernels converged and there 
was no benefit in including additional coefficients in the model. The first-order kernel shown in Figure 7 is 
composed of 257 wavelet functions and has a memory of 2 seconds. The symmetric second-order kernel 
is represented in terms of 561 unique wavelet coefficients and has a memory of 1 second. 

The output predicted by the identified first-order kernel is shown in Figure 8. Also shown is a com- 
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Collective Aileron Input Accelerometer Response 



Time (s) Time (s) 

Figure 5 : Collective aileron input and accelerometer response (right wing forward, inside wing 
fold) at flight condition MachO. 85, 15,000 ft. 


Wavelet-Filtered Data Residual Data 



Time (s) Time (s) 

Figure 6: Morlet-filtered data and residual. 

parison of the Morlet-filtered data and the identified first-order data over a small window in time, corre- 
sponding to the excitation of the 14 Hz mode. The identified output closely matches the Morlet-filtered 
data with some amplitude difference. 

The response predicted by the identified second-order kernel is depicted in Figure 9. The predicted 
second-order response is mostly concentrated in the 6 — 9 second time range. A detailed analysis of the 
residual data revealed significant 12 and 14 Hz responses corresponding to input frequencies of 6 and 7 
Hz, respectively. This occurred in the time range of 6 — 9 seconds and is clearly indicative of a quadratic 
nonlinearity. The fact that this response is evident in other data sets demonstrates that it is indeed due to 
a nonlinearity and not merely a linear response to external disturbances. As shown in the zoomed-in plot 
in Figure 9, the second-order kernel is able to predict this nonlinear response to a reasonable degree of 
accuracy. 
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Figure 7: Identified first and second-order Volterra kernels. 


Predicted Linear Response 
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Figure 8: Linear output predicted by the identified first-order kernel. 

CONCLUSIONS 



Predicted Linear Response 



Data decompositions are investigated for more reliable stability estimates while accounting for uncer- 
tainty, noise, and nonlinearity. Error bounds on the estimates are found to improve with nonstationary data 
processing methods which extract localized time and frequency content. Adaptive data decompositions 
offer a choice of basis functions for enhanced information extraction. 

Morlet wavelet filtering provides a convenient decomposition of the data into a linear component and 
a residual composed of nonlinear response and noise. Volterra operators are designed with orthonor- 
mal, piecewise-polynomial multiwavelets to approximate the first and second-order kernels of nonlinear 
aeroelastic systems. This class of multiwavelets is well-suited for the approximation of Volterra kernels 
because the functions are orthogonal and are easily adapted to the finite domains over which the kernels 
are supported. 
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Predicted Nonlinear Response Predicted Nonlinear Response 



Figure 9: Nonlinear output predicted by the identified second-order kernel. 


A first-order Volterra kernel is identified from filtered data to model the linear dynamics. The nonlinear 
dynamics are then modeled in terms of the second-order Volterra kernel, identified from the residual data. 
An example was given demonstrating the potential of this approach for obtaining reasonably accurate 
models of nonlinear aeroelastic systems with a relatively low number of wavelet basis functions. 

Dryden Flight Research Center 

National Aeronautics and Space Administration 

Edwards, California, April 14, 2003 
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